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Abstract: The optimization-based quantitative determination of multianalyte concentrations 
from biased biosensor responses is investigated under internal and external diffusion-limited 
conditions. A computational model of a biocatalytic amperometric biosensor utilizing a 
mono-enzyme-catalyzed (nonspecific) competitive conversion of two substrates was used 
to generate pseudo-experimental responses to mixtures of compounds. The influence of 
possible perturbations of the biosensor signal, due to a white noise- and temperature-induced 
trend, on the precision of the concentration determination has been investigated for different 
configurations of the biosensor operation. The optimization method was found to be 
suitable and accurate enough for the quantitative determination of the concentrations of 
the compounds from a given biosensor transient response. The computational experiments 
showed a complex dependence of the precision of the concentration estimation on the 
relative thickness of the outer diffusion layer, as well as on whether the biosensor operates 
under diffusion- or kinetics -limited conditions. When the biosensor response is affected by 
the induced exponential trend, the duration of the biosensor action can be optimized for 
increasing the accuracy of the quantitative analysis. 
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1. Introduction 

Amperometric biosensors are analytical devices that measure the changes in the output current on the 
working electrode, due to the direct oxidation or reduction of the products of biochemical reactions [1,2]. 
The device specificity to particular analytes (substrates) is achieved by applying a biological material, 
usually an enzyme [3,4]. The calculation of analyte concentration from the biosensors signal (response) 
is rather simple in the case of a linear dependence of the signal on the substance concentration (linear 
calibration) [2,4]. The problem becomes more complex in the case of non-linear calibration or in the 
presence of mixtures of substances, producing the biosensors response [5-7]. 

Multianalytes have been successfully analyzed by biosensors using different multivariate approaches 
(e.g., partial least squares regression, principal component analysis) [8-10], artificial neural 
networks [11-15] and optimization approaches [16-18]. In those analytical systems, several enzymes or 
even several enzyme electrodes were used when components of the mixtures had an additive effect on 
biosensor response. 

The multianalytes determination becomes even more complex if the biosensors response is perturbed 
by noise, e.g., white noise, sinusoidal power electrical noise or if the biosensor response is biased, e.g., 
by temperature change [19-21]. While in most electrochemical systems, the electrode noise and signal 
trend elimination are practically impossible, a comprehensive study of the background current provides 
valuable insight for analytical system design [19-21]. 

Recently, the problem of multianalytes determination by reverse problem solving has been explicitly 
formulated and successfully applied to optimize the calculation of multianalyte concentration using 
a response of a biocatalytic amperometric biosensor utilizing a mono-enzyme-catalyzed (nonspecific) 
conversion of multi-substrates [22]. The influence of the white noise-, as well as temperature-induced 
trend on the calculation of the analyte concentration has been also investigated at internal diffusion 
limitations by ignoring the external mass transport by diffusion. However, in practical biosensing 
systems, the diffusion of materials outside the enzyme region is of crucial importance for the biosensor 
response [23,24]. 

Real biosensors contain as an outer membrane a thin layer of porous or perforated polyvinyl 
alcohol, polyurethane, cellulose, latex or other material [2,4]. The outer membrane imposes an 
additional diffusion limitation to accessing the substrate to the catalytic layer and, therefore, increases 
the stability and prolongs the calibration curve of a biosensor. The external mass transport should 
be taken into consideration, even if no outer membrane is used and the Nernst diffusion layer 
approximation is applied [25,26]. In those cases, the mass transport outside the enzyme layer is usually 
modeled by an external diffusion layer [4,24-27]. In this work, the optimization-based method of the 
quantitative analysis of the biosensor response [22] is extended by taking into consideration the external 
mass transport. 
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The aim of this work is to investigate the effect of the external diffusion limitation on the precision 
of the determination of the multianalyte concentrations from the biased response of a biocatalytic 
amperometric biosensor utilizing a mono-enzyme-catalyzed multi- substrate conversion. The influence 
of the white noise-, as well as temperature-induced trend one the precision of multianalyte determination 
is also investigated. The investigation was carried out at very different catalytic conditions, thicknesses 
of the diffusion layer and types of signal noise. 

Since a comprehensive analysis of a chemometric technique to be used for quantitative analysis 
usually requires a lot of input data, the biosensor responses to mixtures of compounds were simulated. 
On the other hand, computer simulation is usually much cheaper and faster than real experiments. 
Pseudo-experimental responses to mixtures of two compounds were numerically simulated using a 
computational model of the amperometric biosensor utilizing a mono-enzyme-catalyzed (nonspecific) 
conversion of two substrates [22]. The task of our investigation was to analyze the effect of 
temperature-induced biosensors response drift to analyte concentration calculation. Therefore, the 
simulated responses of the biosensor were perturbed by temperature-induced drift. The simulated 
responses to different concentrations of the substrates were used to extract the dependence of the transient 
output signal on the substrate concentrations. The resulting information was then used to determine the 
concentrations of the substrates from the testing simulated measurements. This work does not analyze 
the determination of individual analytes in the mixture. Rather, the concentrations of known different 
analytes were determined with the biosensor showing unspecific substrate conversion. 

The numerical simulation is based on a two compartment mathematical model involving coupled 
reaction-diffusion equations [27,28]. The model comprises three regions: an enzyme layer, where 
enzymatic reaction, as well as the mass transport by diffusion takes place, a diffusion limiting region, 
where only the diffusion takes place, and a convective region. The numerical simulation was carried out 
using the finite difference technique [26,28,29]. The computational model of the biosensor was validated 
using known analytical solutions for mono-enzyme single substrate amperometric biosensors [27] and 
against published experimental data [30] . 

The numerical experiments showed that the precision of the concentration estimation significantly 
depends on the relative thickness (Biot number) of the outer diffusion layer. The accuracy of 
concentration estimation also depends on the biosensor response being under either the diffusion or 
the enzyme kinetics control, i.e., on the diffusion module. The obtained dependencies can be applied in 
the development of smart analytical systems providing a high quality quantitative analysis of mixtures. 

2. Mathematical Model 

We consider a mono-enzyme dual biosensor (two-substrates) utilizing the Michaelis-Menten 
kinetics [1,3,30], 

E + SiT ^ESi vE + Pi (la) 

E + S2 ^ ^ ESa ^E + P2 (lb) 

k.2 
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where E denotes the enzyme, Si and S2 are the substrates to be determined, ESi, ES2 stand for the 
enzyme- substrate complexes. Pi, P2 are the reaction products and fc+i, A;_i, A;+2, ^-2, ^3, ^4 are the 
kinetic constants. 

When two substrates (Si and S2) react with a single enzyme, E, without the formation of any two 
substrate complex, and the substrates do not combine directly with each other, then each substrate in a 
mixture of Si and S2 acts as a competitive inhibitor of the others [31]. Particularly, 5a-androstan-3-one 
and 5a-androstane-3,16-dione are a special case of the substrates. Si and S2, and 20/3-hydroxy 
steroid-NADoxidoreductase (EC 1.1. 1.53) is a sample of the enzyme, E [30]. 

The reactions in the network given by Equation (1) are usually of different rates [1,31]. To sidestep 
the problems arising due to the large difference of timescales, the quasi-steady-state approach (QSSA) 
is often applied [30]. According to the QSSA, the concentration of the intermediate complex does not 
change with time; the reaction network (Equation (1)) reduces. 



Si 



E 

S2 ^Pi 



(2) 



The amperometric biosensor is treated as an electrode, and a relatively thin layer of an enzyme 
(enzyme membrane) is applied onto the electrode surface. The biosensor model involves three regions: 
the enzyme layer, where the biochemical reactions (Equation (1)), as well as the mass transport by 
diffusion takes place, the diffusion layer, where only the mass transport by diffusion of the substrates, 
as well as the products takes place, and a convective region, where the concentration of the substrates 
and products remains constant. Figure 1 shows the principal structure of the biosensor, where a; = 0 
represents the electrode surface and x = d corresponds to the boundary between the enzyme layer and 
the bulk solution; and, 6 is the thickness of the external diffusion layer. 

Figure 1. Principal structure of the biosensor. 
X 



Bulk solution 



Diffusion layer 
(diffusion) 



Enzyme membrane 
(reaction and diffusion) 



T 

M 
1\ 



Electrode 

Assuming a symmetrical geometry of the electrode and a homogeneous distribution of the 
immobilized enzyme in the enzyme layer of a uniform thickness, the mathematical model of the 
biosensor action can be defined in a one-dimensional-in- space domain [26,27]. 



2.1. Governing Equations 



Coupling the enzyme-catalyzed reactions given by Equation (1) in the enzyme layer with the 
one-dimensional-in-space diffusion leads to a system of the reaction-diffusion equations. 
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dt 
OR, 



D 



D 



dt dx^ ' 1 + Si,e/Ki + S2,e/K2 

where x and t stand for space and time, respectively, Si^eix, t) and Pi^eix, t) are the molar concentrations 
of the substrate, Sj, and the product, P^, in the enzyme layer, respectively, Vi is the maximal enzymatic 
rate, Ki is the Michaelis constant, d is the thickness of enzyme layer, -Ds^.e and Dp-^e are the diffusion 
coefficients, Vi = k^Eo, V2 = k^Eo, Ki = (/c_i + ks)/k+i, K2 = {k-2 + k4)/k+2, and Eq is the total 
concentration of the enzyme, i = 1,2. 

Outside the enzyme layer, only the mass transport by diffusion of the substrates and products 
takes place. 



+ 



1 + SijK, + S2,e/K2 
{V,/K,)S,,e 



(3) 



1,2 0<x<d, t>0 



dt 



D 



s„b- 



D 



dx^ 



(4) 



dt ^^"^'dx^ 

where 5 is the thickness of the external diffusion layer, Si^^x, t) and Pifi{x, t) are the concentrations of 
the substrate, Sj, and the product, Pj, in the bulk and Ds^^h and Dp- }, are the diffusion coefficients of the 
species in the bulk solution, i = 1,2. 



i,b 



1,2 d < X < d + 5 t>0 



2.2. Initial and Boundary Conditions 

The biosensor operation starts when both substrates (Si and S2) appear in the bulk solution {t = 0), 

Si^e{x, 0) = 0, Pi^e{x, 0) = 0, 0<x <d 

5i,fe(x,0) = 0, Pi,fe(x,0) = 0, d<x<d + 6 (5) 
Si,b{d + 6,0) = Si,o, Pi,b{d + 6, 0) = 0, z = l,2 

where Si^ is the concentration of the substrate, Sj, in the bulk solution, i = 1,2. 

At the electrode surface (x = 0), due to the electrode polarization, the concentrations of the reaction 

products (Pi and P2) are permanently reduced to zero, while for the non-ionized substrates, their fluxes 

are assumed to be zero (t > 0) [27], 

dS^ 



P,,e(0,t) = 0, Ds.,e- 



dx 



0, 



1,2 



(6) 



x=0 



Assuming a well-stirred buffer solution leads to the constant thickness of the diffusion layer, as well 
as the constant concentration above that layer during the biosensor operation {t > 0), 



S,,b{d + 6,t) = S,,o, P^,b{d + 6,t) = 0, 



1,2 



(7) 



On the boundary between two adjacent layers, the merge conditions are defined for the substrates, as 
well as the products (t > 0), 



D 



dS, 



Si,e 



D 



Pi,e- 



dx 
dP^ 

dx 



D 



dS,, 



x=d 



x=d 



S^,b 



D 



P,.,b- 



dx 
dl\b 
dx 



x=d 



x=d 



Si^e{d,t) = Si^b{d,t) 
P,,e{d,t) = P,,b{d,t) 



(8) 



Sensors 2014, 14 



4639 



The diffusion layer {d < x < d + 6) may be treated as the Nernst diffusion layer [28]. According to 
the Nernst approach, the diffusion layer of the thickness, 5, remains unchanged with time. 



2.3. Biosensor Response 

The measured anodic or cathodic current is usually assumed as the response of an amperometric 
biosensor. When modeling amperometric biosensors, due to the direct proportionality of the current to 
the area of the electrode surface, the current is often normalized with that area [26,27]. The density, J(t), 
of the biosensor current at time t depends upon the fluxes of the both products at the electrode surface 
and can be expressed explicitly from the Faraday and the Fick laws [31], 

, /oo = lim I{t) (9) 

where rii and n2 are the numbers of electrons involved in charge transfers in the corresponding 
electrochemical reactions at the electrode surface, I^o is the density of the steady-state biosensor current 
and F is the Faraday constant, F = 96, 486 C/mol. 



I{t) = niFDp,^ 



dx 



+ n2FDp^^ 



OR 



2,e 



.x=n 



dx 



2.4. Characteristics of Biosensor Action 

The diffusion module or Damkohler number essentially compares the rate of the enzyme reaction 
(Vi/Ki) with the mass transport through the enzyme-loaded layer {Dsi,e/d'^) [27,32], 

9 Vid'^ 

$? = -^, 1 = 1,2 (10) 

where $f is the dimensionless diffusion module corresponding to the i—th reaction of the network given 
by Equation (1). 

If the diffusion module is less then unity, then enzyme kinetics (reaction rate) controls the biosensor 
response. The response is controlled or limited by diffusion when the module is greater than unity [27]. 
At intermediate values of the diffusion module, the biosensor operation is of a mixed control. The 
quantitative analysis of mixtures may notably depend on whether the biosensor response is under 
diffusion or enzyme kinetics control [12,13,22]. 

The Biot number is another dimensionless parameter widely used to indicate the internal mass transfer 
resistance to the external one [4,33], 

where (3i and (32 are two Biot numbers corresponding to the diffusion of two substrates. 



2.5. Signal Trend and Noise 

In real-life experiments, the reaction rate can be affected by the temperature change that induces 
a trend of the output signal of the sensor [19,34]. The trend factor, R{t), is often expressed by the 
Arrhenius equation, which defines the dependency of the reaction rate on the temperature [35]. The 
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trended signal, Irif), can be expressed as the multiplication of the initial output, I{t), by a trend factor, 

R{t),t > 0, 

Irit) = I{t) X R{t) (12) 

Assuming that the temperature, T, changes linearly with time (T = Tq + at) and -R(O) equals unity, 
the trend factor, R{t), is expressed as follows: 

Kit) = exp f ^ X (13) 

where Ea is the activation energy, R = 1.98 (cal/(K ■ mol)) is the gas constant, T is absolute temperature. 
To = 298 K, and a (K/s) is the coefficient of proportionality [22,35]. The zero value of a corresponds to 
a signal effected with no trend, Irit) = lit) at a = 0 and t > 0. 

Besides the signal trend, the biosensor response can be also affected by an unpredictable noise [20,2 1 ] . 
The measurements, I Nit), of the noisy signal can be modeled by adding a white Gaussian noise to the 
noise free biosensor response, lit), 

/^(t)=/(t)+e (14) 

where ^ ~ A/'(0, crlit)) stands for the random number, generated following the Gaussian distribution 
with zero mean and standard deviation equal to (j/(t) [22]. Then, the trended noisy output signal IxNit) 
can be modeled as follows: 

lTNit) = Iit)Rit)+^ (15) 

The noise can be signal-dependent or signal-independent [20]. The electrical background noise is 
one of the major limitations of the amperometric detection of the smallest possible concentrations of 
electrochemically active materials [36]. The background noise is usually signal independent. This 
work concentrates on the determination of relatively high concentrations of substrates in comparison 
with the corresponding Michaelis constant. The multiplicative signal-dependent noise is used in the 
modeling, since the measurement errors in such applications are usually proportional to the measured 
signal value [20,36]. The influence of the additive noise can be also investigated; e.g., the influence of 
both kinds of noise was compared in [22] for a similar model, but without an external diffusion layer. 

2.6. Numerical Simulation 

Due to the nonlinearity of the governing Equation (3), the initial boundary value problem. 
Equations (3)-(8), can be analytically solved only for a very specific set of the model parameters [27,28]. 
Because of this, the problem. Equations (3)-(8), was solved numerically by applying the finite difference 
technique [26,28,29]. In the space coordinate, x, the enzyme and the external diffusion layers were 
divided into the same number of subintervals of equal length. A uniform discrete grid was also introduced 
in the time coordinate, t. An explicit finite difference scheme has been built as a result of the difference 
approximation of Equations (3)-(8) [22,29]. The numerical simulator of the biosensor action has been 
programmed in the C-i-i- language [37]. 

Explicit difference schemes have a convenient algorithm of the calculation and are simple for 
programming [26,28]. The explicit scheme used for Equations (3)-(8) was conditionally stable 
and converges with the rate 0(r + max (/ii, /i2)^), where hi = d/N, /12 = 6/N. In the 
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simulation, both layers were discretized with the same number = 200 of grid points. was 
constant in the simulation of all the responses, while the time step size, r, was recalculated for each 
simulation. To make the difference scheme stable, the step size, r, was found from the sufficient 
stability conditions, r max g, Dp. e) < h\/A, r max (1^5. e, -Dp^^e) < ^i/^j i = 1,2, and 
r max {Vi, V2) < {K^ + K2)/2 [26,38]. 

The mathematical, as well as the corresponding computational models of the biosensor were validated 
using known analytical solutions for mono-enzyme single substrate amperometric biosensors [27]. When 
the concentration of the substrates is low in comparison with the corresponding Michaelis constant 
{Si^e < Sifi <ti Ki), the nonlinear reaction terms in Equation (3) simplify to those of the first order, 
{Vi/Ki)S,,i = l,2. 

For validating the model in the opposite case of the substrate concentrations, a high concentration, 
Si,o, of the substrate, Sj, was used together with a zero concentration of another substrate, Si 0 ^ Ki, 
5*2,0 = 0 and 6*2,0 ^ K2, S'1,0 = 0. In both these cases of the extreme concentrations of the substrates, 
the initial boundary value problem. Equations (3)-(8), reduces to a linear one and can be solved 
analytically [22,26,27]. Reducing the thickness of the diffusion layer {5 — t- 0), the two compartment 
model. Equations (3)-(8), approaches a known one compartment model [22], the numerical solution of 
which was also used for validating the numerical simulation of the problem. Equations (3)-(8). 

The calculated values of the reaction term Equation (3), as of the total rate of Reaction (1) 
at steady-state conditions was also validated using published experimental data [30], accepting 
5a-androstan-3-one and 5a-androstane-3,16-dione as the substrates. Si and S2, and 20/3-hydroxy 
steroid-NAD oxidoreductase (EC 1.1.1.53) as the enzyme, E. 

The following values of the model parameters were constant in all the numerical experiments: 

D5.,e = /^p„e = 3 X 10-^ cmVs, = ^p.fc = 6 X 10-6 cmVs ^^^^ 

Ki = 10'^ M, rf = 0.02cm, rn = 1, i = 1,2 
The membranes were characterized by a diffusion coefficient of substrates. For the calculation, the 
diffusion coefficient, SxlO^^ cm^/s, was used, which is typical for the enzyme membranes ofbiosensors 
and is two times less than that in the buffer solution [39]. 

In order to investigate the impact of the external diffusion limitation on the precision of the 
determination of the multianalyte concentrations and to avoid results dependent on the predefined values 
Equation (16), the biosensor response was simulated at different values of the thickness, 5, of the 
diffusion layer chosen, so that the Biot numbers Equation (11) would vary in a wide range (from 0.4 
up to oo). The infinite Biot number corresponds to the zero thickness of the diffusion layer or to the 
extremely low diffusivity of materials in the diffusion layer. Although the thickness, S, of the diffusion 
layer was varied in the investigation, similar responses could be also simulated by changing the diffusion 
coefficients, Dsi^ and -Dp.,;, {i = 1, 2), of materials in the external diffusion layer, keeping the same Biot 
numbers. To minimize the results dependent on the catalytic properties of the biosensor, the response 
was simulated at different value of diffusion modules Equation (10) by changing the maximal enzymatic 
rates, Vi and V2, so that the response would be in both enzyme kinetics and diffusion limitations. 

The numerical solutions of the model. Equations (3)-(8), were compared with the analytical solutions 
at different values of Vi and Si^, i = 1,2. The relative difference between the numerical and analytical 
solutions was less than 1% when applying the stability conditions defined above. 
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The rate of the trend is characterized by the activation energy, Ea, and the coefficient, a, of 
proportionality. The simulated biosensor responses were affected by different trends using different 
values of the latter parameters. Values of the trend parameters used in the calculations are presented in 
Table 1 , where the first column stands for the denotation of the trend, the second one for the value of 
activation energy, Ea and the last column for the value of the coefficient, a, of proportionality. These 
values were selected to be the same as in [22] and to model trends similar to the temperature-induced 
experimental ones [40] . 

Table 1. The values of the parameters of the exponential trend Equation (13) substantiated 
in [22]. 



Title 


Ea (cal/mol) 


a(K/s) 


Trend 1 


6,000 


3.33 X 10-3 


Trend 2 


24,000 


3.33 X 10-3 


Trend 3 


24,000 


6.66 X 10-3 



The evolution of the typical biosensor responses is illustrated in Figure 2, where the left plot shows 
the signals affected by the exponential trend of a different rate and the right plot the response affected by 
the trend, entitled Trend 2, in conjugation with the white Gaussian noise at cr = 0.05. The simulation was 
performed at the maximal enzymatic rate Vi = 0.5 |J.M/s for the first substrate and V2 = 10 Vi = 5 |J.M/s 
for the second one; the concentration of the first substrate was Si^ = 3.2Ki = 0.32 mM and 
5'2,o = 4Si,o= 12.8/^2 = 1.28 mM of the second substrate, and the thickness of the outer diffusion layer 
was (5 = 1 mm. Values of the other parameters were as defined in Equation (16). At the values used in the 
simulation, the Biot numbers were /3i = P2 = 0.4 and the diffusion modules: ^ 0.67, <I>1 = 10^1^ 6.7. 

Figure 2. The dynamics of the biosensor current affected by the exponential trends (Left) 
and in conjugation with the white Gaussian noise (Right, a = 0.05) simulated at the 
following values of the model parameters: Vi = 0.5 [xM/s, V2 = 5 |J.M/s, 5*1,0 = 0.32 mM, 
5*2,0 = 1-28 mM. The other parameters are as defined in Equation (16). 



n No Trend 

Trend 1 




0 600 1200 1800 2400 3000 3600 0 600 1200 1800 2400 3000 3600 

t,s t,s 



3. Optimization Problem 

We consider the quantification of mixtures of two substrates from the response of the mono-enzyme 
biosensor. Similar problems for simpler models have been already considered in [18,22]. In this work. 
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we consider the quantitative determination of two substrates acting as competitive inhibitors of each 
other at the external diffusion limitation. As in [22], the measurements are supposed to be corrupted 
by a white Gaussian noise, and the reaction rate is to be affected by the temperature variance [7,19]. 
Pseudo-experimental simulated data were used for investigating the effect of diffusion limitations on the 
the precision of the concentration estimation. A similar approach to using simulated biosensor responses 
has been applied for training a neural network to perform quantitative analysis of the mixtures [12,13]. 

3.1. Statement of the Optimization Problem 

Let W = {wi, . . . , Wn) be a sequence of test measurements of the biosensor current at discrete time 
moments, ti, . . . , tn, and let Z{c) = {z{ti, c), . . . , z{tn, c)) be a sequence of the corresponding values 
of the response numerically simulated using the model. Equations (3)-(8), where c = (ci, C2) stands for 
the concentrations of the substrates. Si and S2, which are the subject to be evaluated, and z{ti, c) is the 
simulated biosensor current at time moment ti, i = 1, . . . ,n. 

The concentrations c = (01,02) can be evaluated by tuning the values of the concentrations with 
respect to the conformity of theoretical measurements, Z(o), to the given sequence, W, of the measured 
biosensor currents. In order to perform the tuning, the least squares approach has been applied in [18], 
where the relevant problem has been described as: 

n 

0 = argminy^iziti, 0) - WiY (17) 

i=l 

where C denotes the feasible region of the values of concentrations c = (ci, C2). 

When the biosensor output is affected by the exponential trend Equation (13), the mathematical 
optimization problem given by Equation (17) should be reformulated as: 

n 

(c, a) = arg min ^{z(ti,o) R{ti) - Wif , (18) 

{c£C, a_<a<a+} ' ' 
1=1 

where the decision parameter, a, is supposed to be unknown in Equation (13) and the subject is to be 
evaluated within the feasible interval [a_, a+]. 

The optimization problem defined in Equation (18) is a generalization of the optimization problem 
given by Equation (17), and both of them coincide when R{t) = 1 or; in other words, when Ea = 0 
and/or a = 0 in Equation (13). 

3.2. Solution of the Optimization Problem 

The optimization problems stated above are difficult to solve, due to mathematically unprovable 
convexity and uni-modality properties of the objective function, thus limiting the applicability of local 
optimization algorithms, and are highly computational resource-consuming to find a single value of 
the objective function. The approach to apply a surrogate objective function has been proposed in [22], 
where the surrogate objective function, c), is defined as bi-linear interpolant of the function, z{ti, c), 
assuming that the theoretical biosensor responses, z{ti, c), are already simulated using a discrete set, C", 
of concentrations c = (ci, C2) for alH = 1, . . . , n. 
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In the optimization problem given by Equation (18), the objective is the squared error of the measured 
biosensor response {wi , . . . , m„) at discrete time moments 1 1 , . . . , t„ with respect to known (simulated or 
measured experimentally) biosensor responses Z{c) obtained at the feasible region, C, of the substrate 
concentrations. The decision variables are the concentrations, Ci and C2, to be determined, and the 
feasible region, C, is defined by their intervals. For the optimization, a state-of-the-art algorithm was 
used [22]. The maximal time, t„, is chosen so large that all the biosensor responses, Z{c), approach the 
steady state theretofore. 

The feasible region for the concentrations c = (ci,C2) has been chosen to be C = [3.2, 12.8]^, 
assuming that the concentrations are dimensionless and normalized with respect to the corresponding 
Michaehs constant as follows: Cj = Si^/Ki, i = 1, 2. The biosensors cahbration curve is usually linear 
at the substrate concentrations lower than the Michaelis constant, and then, the determination of the 
analyte concentration from the biosensor response is rather simple [1,2,4]. Because of this, the feasible 
region was chosen so that Sifi > Ki, i = 1,2, when the concentration determination is more complex. 

The entire domain, C, has been discretized with the uniform grid with a step size equal to 0.1, the 
nodes of which (97^ = 9, 409 in total) comprise the discrete set, C. The theoretical biosensor response 
values, z{ti,c), have been computed for every combination, (01,02) G C, assuming that to = 0; 
tj+i — ti = Is; i = 1, . . . ,n — 1; n = 3,600. Such a predefined duration (t„ = 3,600s) of the 
biosensor operation is sufficient to reach the steady-state current using the thickest diffusion layer and 
lowest reaction rate that was used in the computational experiments. 

The computation of a pseudo-experimental response Z(o) = {z{ti, 0), . . . , z{tn, 0)) for a particular 
pair c = (ci,C2) of the concentrations requires from around 15 to around 35 min (depending on the 
thickness of the diffusion layer and the reaction rate) using one core of an Intel Xeon X7550 (2 GHz) 
processor. Therefore, the computation of the responses for the whole set, C , of the concentrations is time 
consuming or even impossible within reasonable time. On the other hand, the simulation of different 
biosensor responses can be considered as independent tasks and performed in parallel by assigning 
them to different computing resources. For this purpose, the framework for computational simulation 
of a large set of biosensor responses using high performance computing facilities has been developed. 
The framework is based on the master-slave strategy, where one of the processing units (the master) is 
considered to be responsible for the management of the queue of simulations assigned to be performed, 
and the others (the slaves) performs computational work assigned by the master [41]. The framework has 
been used on a high performance facility, consisting of 32 Intel Xeon X7550 processors with eight cores 
each (256 cores in total). Due to the homogeneity of the processing units and the relatively small costs 
required for communication between them, the speedup of the computations was almost linear (almost 
equal to the number of processing units that were used), thus allowing one to simulate the responses, 
Z{o), for a whole set, C , of pairs c = (ci, 02) of the concentrations within reasonable time. 

Since the simulation of a single biosensor response requires from 15 up to 35 min when using a single 
core of an Intel Xeon X7550 processor, the simulation of 9, 409 such responses would require more 
than 2, 352 h, whereas the latter computations using 256 cores of Intel Xeon X7550 processor require 
less than 10 h. This indicates the almost linear speed-up of the computations or, in the other words, a 
reduction of the runtime by almost 256 times. 
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The obtained pseudo-experimental responses, Z{c) (c G C), have been used to derive the values of 
the surrogate objective function, c), by applying the bi-linear interpolation approach. Because of the 
availability of the analytical expression of the gradient of z(t, c), a standard gradient descend method 
can be applied to minimize z{t,c). 

The optimization-based method for the evaluation of the concentrations has been implemented in 
a MATLAB environment using the Multi-Startstrategy with an efficient local minimization function, 
fmincon, which is included into the optimization toolbox as a subroutine [42]. A user-defined formula 
was used for the gradients of the objective functions. Equations (17) and (18). The local search has been 
performed three times using different initial solutions, randomly generated within the feasible region. 

The simulation of the biosensor responses is much more time-consuming than the optimization 
process; therefore the response simulator has been programmed in C-i-i-, including distributed memory 
parallel programming libraries. 

The standard stopping criteria for the local search has been defined by a tolerance TolFun = 10~^ 
of the function value [42]. Since the scale of the biosensor response varies depending on the model 
parameters, the same function value corresponds to the different discrepancy between the simulated and 
the given measurements of the biosensor response. To unify the scale of the biosensor response and the 
scale of the objective function values, the measurements of the biosensor response have been multiplied 
by a scalar, chosen with respect to making the maximal objective function value equal to 10^. 

4. Results and Discussion 

The optimization-based method for the evaluation of the concentrations (see Section 3.2) has been 
used to investigate the precision of the evaluation of the concentrations, when the biosensor response 
follows the model. Equations (3)-(9). Each component of the mixture was characterized by the 
individual maximal enzymatic rate differing in an order of magnitude. Without reducing the generality, 
it was assumed that the maximal enzymatic rate, V2, for the second substrate is greater than the rate, Vi, 
for the first one, V2 = lOVi and $2 = 10$i. Similar approaches have been already used in a quantitative 
analysis of mixtures in [12,13,22]. 

The quality of the quantitative analysis is influenced by whether the biosensor response is under the 
diffusion or the enzyme kinetics control, and the concentration estimation is usually more accurate for a 
substrate corresponding to a greater diffusion module than for another substrate corresponding to a lower 
diffusion module [12,13,22]. Because of this, the maximal enzymatic rates, Vi and V2, were chosen, so 
that the biosensor response can be controlled by the enzyme kinetics for the first component ($^ < 1) 
and by the mass transport for the second one ($2 > 1), Vi = 0.5 |J.M/s, V2 = 10 Vi = 5 |J.M/s, $1 ^ 0.67, 
$2 = 10$^ ^ 6.7. 

The precision of the evaluation has been determined by attempting to evaluate the set of 64 pairs 
c = (ci,C2) of predefined concentrations for which the responses were simulated (see Figure 3) and 
calculating the average precision of the estimates c = (01,02) obtained by the optimization-based 
method. The precision of a single estimate of the concentration has been determined by the 
relative error. 
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Figure 3. The set of the dimensionless concentrations, Ci and C2, of the substrates. Si and S2, 
used in the investigation, as well as the noise-free biosensor responses to these concentrations 
simulated at and Vi = 0.5 ]xM/s, ¥2 = 5 \iM/s and (3 = 0.4. 
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The sensitivity of the estimates of the concentrations is illustrated by Figure 4, where contour lines 
of the objective function given by Equation (17) for the noise-free biosensor responses were simulated 
at Vi = 0.5 ia.M/s, V2 = 5 |J.M/s assuming zero thickness {(3 = 00) of the external diffusion layer, as well 
as the relatively thick (/3 = 0.4) diffusion layer, where the true dimensionless concentrations (ci, C2 ) of 
the substrates were equal to (8, 8). As one can see in Figure 4, the objective function for both values of 
/3-parameter is unimodal in the region of interest, and the minimum point lies in a wide flat valley, which 
is curved and can negatively influence the convergence of the gradient-based algorithm when (3 = OA. 

Figure 4. Contour lines of the objective function given by Equation (17) for the noise-free 
biosensor responses simulated at Vi = 0.5 \^M/s, V2 = 5 \lM/s assuming zero (/3 = 00) and 
a relatively thick (/3 = 0.4) thickness of the external diffusion layer, (ci, C2) = (8, 8). 
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The maximal time t„ = 3,200 s until which the responses were measured was chosen so large 
that all the biosensor responses approach the steady state. Probably, a reasonable quality of the 
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concentration determination can be achieved from shorter responses when the steady state is not reached. 
In order to investigate the influence of the duration of the biosensor operation on the precision of 
the concentration evaluation, the concentration evaluation has been performed by analyzing biosensor 
responses of different duration of the biosensor operations;, G {200,400,800,1,600,3,200}, 
first seconds of the operation, < tn- The duration tm = tn = 3,200 s corresponds to the 
steady-state response. 

In order to investigate the impact of the thickness, S, of the external diffusion layer, the biosensor 
responses were simulated at different values of S chosen so that the Biot number (3 {(3 = (3i = (32) would 
vary in a wide range (from 0.4 up to oo), assuming the thickness, d, of the enzyme layer to be constant, 
as defined in Equation (16). The infinite Biot number (/3 = oo) corresponds to the zero thickness of the 
diffusion layer. 

Since the concentrations from the noise-free biosensor responses can be evaluated very 
precisely [13,18] and the measurements in real experiments are often affected by an unpredictable 
noise [19-21], only the noisy responses modeled by adding either an entirely white Gaussian noise 
Equation (14) or that noise in conjugation with the induced trend Equation (15) were used for 
investigating the precision of the concentration evaluation. 

4.1. Impact of Biosensor Response Time 

The response time is an important characteristic of biosensors. In various applications, it is important 
to have the response time as short as possible [1,3]. The response time of enzymatic biosensors depends 
to a great extent on the diffusion processes [4,26]. Particularly, an increase in the thickness of the 
external diffusion layer prolongs the steady-state response time. When the steady-state response time is 
relative large the biosensor current measured much earlier than the steady-state achieved can be used as 
the biosensor response. Having only partial data on the evolution of the biosensor current can notably 
decrease the precision of the quantitative analysis of the analyte [22]. Because of this, it is important 
to know how short the duration of the biosensor operation can be for the substrate concentrations to be 
suitably evaluated. 

Computational results showed that for all practical values of the thickness, 5, of the outer diffusion 
layer, the concentrations can be evaluated with a relative error lower than 0.05 from the measurements of 
the longest biosensor response {tm = 3, 200 s) corrupted by the white Gaussian noise with a = 0.05 (see 
Equation (14)). However, as one can see in Figure 5, the precision can be impaired using measurements 
of a shorter biosensor response. In Figure 5, different curves stand for the different thicknesses of 
the diffusion layer, corresponding to the different Biot numbers. The impairment of the precision, 
when shortening from t„ to of the biosensor response, has been evaluated by the change of the 
relative error: 

Ae,{tm) = et'^ - Ei (20) 

where €i is the relative error of the estimate of the concentration, q, obtained using measurement of 
the first tn seconds of the biosensor response (see Equation (19)), and e^"^^ is the relative error obtained 
using measurements of the first tm seconds of the biosensor response, i = 1,2. The reference value 
tn = 3, 200 corresponds to the steady-state response. 
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Figure 5. The impairment, A^j, of the relative error, Si, of the evaluation of the 
concentration, q, of the substrate, Sj, versus the duration, tm, of the biosensor operation at 
four values of the Biot number, /3, when the biosensor response is affected by white Gaussian 
noise with (j = 0.05,z = l,2. The other parameters are as in Figure 2. 
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Figure 5 shows that shortening of the biosensor response time, t^, from 3,200 down to 800 s, the 
relative error (A^j, as well as et) of the evaluation of the substrate concentrations arises insignificantly 
or even remains practically the same (e.g., when (3 = oo). At relatively short-term biosensor responses 
(tm = 200 and tm = 400 s), the impairment, Ae^, of the relative error notably increases with decreasing 
the Biot number, (3 (increasing the thickness 6 of the diffusion layer). The impact of the Biot number on 
the accuracy of the concentration evaluation is discussed in detail below. 

Since the maximal enzymatic rate of the second substrate was notably (ten-fold) greater than of 
the first one, the second substrate has more impact on the biosensor response [22,26]. Moreover, the 
concentration, C2, of the second substrate was evaluated notably more precisely than that of the first one 
in all the numerical experiments that were performed. Taking these considerations into account, bellow, 
the investigation is focused on the evaluation only of the concentration, ci, as less precisely determined. 

4.2. Impact of Exponential Signal Trend 

Similar computational experiments were applied for investigating the impact of thickness 5 of the 
outer diffusion layer in conjugation with the temperature-induced trend. The results showed that the 
concentrations of both competitive substrates can be precisely evaluated practically independent of the 
thickness of the diffusion layer, the duration of the biosensor operation and the rate of the exponential 
trend, when the noise-free biosensor responses are analyzed. 

The maximal relative errors, 6i and 62, of the concentrations, ci and C2, were less than 0.02 and 
0.005, respectively. The most intractable situations occurred when the diffusion layer was the thickest 
{(3 = 0.4), the biosensor response was shortest (tm = 200 s) and the signal was affected by the 
exponential trend of the highest rate (Trend 3). 

Figure 6 presents the relative error, si, of the estimate, ci, of the substrate. Si, versus the duration, tm, 
of the biosensor action. The error, ei, was calculated for different values of the Biot number, (3, when the 
biosensor response is affected by white Gaussian noise with a = 0.05 (see Equation (14)) in conjugation 
with the exponential trend (see Equations (13) and (18)). The biosensor current was calculated as defined 
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in Equation (15). In Figure 6, different curves correspond to the different rates of the exponential trend 
(see Table 1), and different plots to the different Biot number, (3. 

Figure 6. The relative error, ei, of the evaluation of the concentration, ci, versus the duration, 
tm, of the biosensor operation at different values of the Biot number, (3, when the response is 
affected by white Gaussian noise with a = 0.05 and the exponential trend of a different rate. 
The other parameters are as in Figure 2. 
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One can see in Figure 6 that the duration of the biosensor transient response noticeably influences the 
precision of the concentration evaluation. The influence depends on the rate of the trend: a higher rate 
leads to a higher impairment of the evaluation precision. On the other hand, the usage of a relatively 
short biosensor response also leads to a significant impairment of the evaluation precision. Therefore, 
the duration of the biosensor action can be optimized when the biosensor response is affected by the 
exponential trend. As one can see from Figure 6, the optimal duration of the response depends on 
the relative thickness (Biot number (3) of the outer diffusion layer: a thicker layer (smaller (3) requires 
a longer biosensor response; it is worth using tm = 400 when (3 = oo, tm = 800 when (3 = 2.5, 
tm = 1, 600 when (3 = 1 and = 3, 200 s (steady-state response), when /3 = 0.4. 

Figure 7 shows the dependence of the precision of the concentration evaluation on the Biot number, (3. 
Different curves corresponds to a different rate of exponential trend, and different images to a different 
duration of the biosensor response, tm G {200, 400, 800, 3, 200} s. The duration of 3,200 s corresponds 
to the steady-state response in all the numerical simulations of the biosensor response. 

One can see in Figure 7 that the relative error, ^i, of the evaluation of the substrate concentration as a 
function of the Biot number can be of a different monotonicity. The error, ei, is a monotonous decreasing 
function of the Biot number, (3, in the case of short-term responses (tm G {200, 400} s), while in the case 



Sensors 2014, 14 



4650 



of long-term responses (t^ G {800, 3, 200} s), the monotonicity of the function, depends on the 

trend rate. 



Figure 7. The relative error, ei, of the evaluation of the concentration, ci, of the substrate, 
Si, versus the Biot number, (3, and the duration, tm, of the biosensor operation when the 
response is affected by white Gaussian noise with a = 0.05 and the exponential trend of a 
different rate. The other parameters are as in Figure 2. 
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Only the dynamics of the biosensor current up to the steady state contains the full information on 
the dynamics of the biosensors operation. A shorter evolution of the biosensor current contains limited 
information used in the evaluation of the substrate concentration. The duration of the biosensor operation 
is especially important at the external diffusion limitation, i.e., a relatively thick outer diffusion layer 
(small values of the Biot number). On the other hand, a longer response is more distorted by the 
exponential trend than a shorter one (see Figure 2). Figure 7 shows that the greatest values of the error. 
El, of the concentration evaluation appear in the case of a thick diffusion layer (small values of /3) when 
short-term responses are used in the evaluation. 

Increasing the thickness of the external diffusion layer creates an additional diffusion limitation to 
the substrates, i.e., leads to a lowering of the substrate concentration in the enzyme layer and, thereby, 
increasing the biosensor sensitivity, as well as prolonging the calibration curve of the biosensor [23,24]. 
This feature can be also noticed in Figure 7. When the response is not affected by the trend and the 
long-term response is analyzed (tm G {800, 3, 200} s), the error, ^i, of the concentration evaluation 
decreases with decreasing the Biot number, (3 (increasing the thickness S of the diffusion layer). 

Figure 7 also clearly shows that at different conditions, the concentrations are more accurately 
evaluated when the response is not affected by a trend rather than affected by a trend. 
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4.3. Impact of Maximal Enzymatic Rates 

In the numerical experiments discussed above, the values Vi = 5 x 10"'' and V2 = 5 x 10"^M/s of 
the maximal enzymatic rates were chosen so that the biosensor response would be under mixed control, 
controlled by the enzyme kinetics for the first component ^ 0.67 < 1) and by the mass transport 
for the second one ($2 ~ 6.7 > 1). In order to investigate the impact of the maximal enzymatic rates 
to the precision of the evaluation, the following two additional configurations of the parameters have 
been investigated: 

• the biosensor response to both substrates is controlled by the enzyme kinetics {Vi = 5 x 10"^, 
1^2 = 5 X 10-^M/s, $2 ^ 0.067, $2 ^ 0.67); 

• the biosensor response to both substrates is controlled by the mass transport (Vi = 5 x 10^^, 
1^2 = 5 X 10-^M/s, $2 ^ 6.7, $2 ^ 67). 

The calculations showed that the difference between the evaluation precisions of different substrates 
changes with the changing of the maximal enzymatic rates, while the concentration, C2, of the second 
substrate, S2, was evaluated notably more precisely in all the numerical experiments discussed above. 

Figure 8 shows the relative errors, ei and 62, of the evaluation of the concentrations of both substrates. 
Si and S2, versus the duration, t^, of the biosensor operation at different values of the maximal enzymatic 
rates, Vi and V2, and a moderate value of the Biot number (3 = 1. The response was affected by white 
Gaussian noise with a = 0.05 and by the exponential trend of a different rate. 

As one can see in Figure 8, when the biosensor response to both substrates is controlled by the 
enzyme kinetics (upper plots), the relative error, ei, of the concentration evaluation for the substrate. Si, 
is noticeably greater than the error, £2, for the second substrate, S2. This is especially noticeable for short 
durations of the biosensor operation, tm < 800 s. A similar property of the concentration evaluation was 
noticed also in the case when the response was under mixed control (see Figure 5). However, when the 
biosensor response to both substrates is controlled by the mass transport (lower plots), the evaluation 
errors are rather similar for the both substrates, especially at short durations of the biosensor operation, 
and even the error, £2, becomes slightly greater than the error, ei, when long-term responses are used in 
the concentration evaluation. 

Figure 8 also shows a noticeably greater error of the substrate evaluation when the biosensor response 
to the both substrates is controlled by the mass transport (lower plots). 

It is known that when the response of a catalytic biosensor is considerably controlled by the mass 
transport, the steady-state current practically does not depend on the maximal enzymatic rate (total 
concentration of enzyme) [26,43]. This feature of amperometric, as well as potentiometric biosensors is 
notably affected by the relative thickness (Biot number) of the outer diffusion layer. Therefore, when the 
biosensor response to both substrates of the mixture is controlled by the mass diffusion ($^ > 1, $2 > 1)' 
the impacts of the both substrates on the response become approximately the same, which surely impairs 
the precision of the concentration evaluations. Choosing the individual maximal enzymatic rates of 
compounds most appropriate for high precision evaluation of the concentrations is the subject of our 
future investigation. 
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Figure 8. The relative errors, ei and 62, of the evaluation of the concentrations, ci and C2, 
of both substrates. Si and S2, versus the duration, tm, of the biosensor operation at different 
values of the maximal enzymatic rates, Vi and V2, and the Biot number (3 = 1, when the 
response is affected by white Gaussian noise with a = 0.05 and the exponential trend of a 
different rate. The other parameters are as in Figure 2. 
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The maximal enzymatic rate, Vi, is actually a product of two parameters: the catalytic constant, k^, 
introduced in Equation (la), and the total concentration, Eq, of the enzyme, and correspondingly, V2 is a 
product of Eq and the constant, k^, introduced in Equation (lb), Vi = k^E^, V2 = k^E^ [1,31]. Since, in 
actual applications, it is usually impossible to change the values of the constants, k^ and k^, the maximal 
rates, Vi and V2, as well as the diffusion modules, and might be changed by changing the total 
concentration, Eq, of the enzyme in the enzyme layer. 



5. Conclusions 



The optimization-based method of the quantitative analysis of the biosensor response proposed in [22] 
is suitable for the evaluation of the concentrations of the competitive substrates from the biased response 
of mono-enzyme amperometric biosensors utilizing the Michaelis-Menten kinetics Equation (1) at 
different diffusion limitations. 

The mathematical model. Equations (3)-(9), can be used to simulate pseudo-experimental biosensor 
responses to mixtures of substrates for evaluating the precision of the quantitative analysis, as well as for 
calibrating an analytical system. 
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If the biosensor signal is not affected by the temperature-induced trend, then the substrate 
concentrations are most accurately evaluated from the biosensor transient response recorded up to a 
steady state. Shortening the duration of the biosensor operation reduces the accuracy of the evaluation, 
especially in the case of a relatively thick outer diffusion layer (small values of the Biot number, /3; 
Figure 5). However, if the biosensor response is affected by the induced exponential trend, then the 
preferable duration of the biosensor action depends on the relative thickness of the outer diffusion 
layer (Figures 6-8). 

At different diffusion limitations and durations of the biosensor operation, the concentrations are more 
accurately evaluated when the response is not affected by a trend rather than affected by an induced 
exponential trend (Figures 7 and 8). 

At different durations of the biosensor operation and rates of the induced trend, the concentrations of 
the substrates are more accurately evaluated when the biosensor response to the substrates is controlled 
by the enzyme kinetics rather than the response being controlled by the mass transport (Figure 8). 
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